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Abstract 

Heterogeneity of economic agents is emphasized in a new trend of 
macroeconomics. Accordingly the new emerging discipline requires one 
to replace the production function, one of key ideas in the conventional 
economics, by an alternative which can take an explicit account of dis- 
tribution of firms' production activities. In this paper we propose a new 
idea referred to as production copula; a copula is an analytic means for 
modeling dependence among variables. Such a production copula pre- 
dicts value added yielded by a firm with given capital and labor in a 
probabilistic way. It is thereby in sharp contrast to the production 
function where the productivity of firms is completely deterministic. 
We demonstrate empirical construction of a production copula using 
financial data of listed firms in Japan. Analysis of the data shows that 
there are significant correlations among their capital, labor and value 
added and confirms that the values added are too widely scattered to 
be represented by a production function. We study four models includ- 
ing the variants of Frank, Gumbel, Clayton copulas and find that the 
non-exchangeable Gumbel copula describes the data distribution and 
the correlations very accurately. 
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I. Introduction 



Recently a new approach to macroeconomics has emerged ( Aoki and Yoshikawa 
( [20071 ), pem Gatti eTall ( [2008| ) by taking serious account of heterogeneity 
of economic agents. Various stylized facts have accumulated to support this 
new perspective, including the distribution of firms' size with a power-law tail. 
The traditional concepts in mainstream economics should thus be replaced by 
new ones: a representative agent, by heterogeneous agents; average, by distri- 
bution; deterministic, by stochastic; mechanical point of view, by statistical 
mechanical point of view; and so forth. The principal goal of the new approach 
is to provide microscopic foundations for macroeconomics. Statistical mechan- 
ics in physics is a discipline to understand macroscopic states of matter using 
microscopic information on atoms and molecules. It is a good successful model 
for the endeavor to bridge between microeconomics and macroeconomics. 
In this paper we pay our attention to production activities of firms. The 



production function is one of basic ingredients in microeconomic theory ( Wick- 
steed (1894), Varian (1992])). It specifies output F of a firm for given input 
factors such as labor L and capital K: 



Y = F{L,K), 



(1) 



where it is assumed that each firm produces goods in an optimized manner 
with its own production function. On the other hand, an aggregated produc- 
tion function is also an important tool to measure production activities at 



national level in macroeconomic theory (Sato (1975l). The GDP, given as the 
total sum of value added of producers over nation, is assumed to be a func- 
tion of aggregated capital and labor for instance. An aggregated production 
function is in principle obtained by summing up all the production functions 
of individual producers. However this naive definition is an impractical way 
to construct the GDP; actually it is estimated through the SNA. There are 
also methodological problems associated with the aggregation, including how 
to define aggregated capital and labor. Thus the production function at a 
macroscopic level has no sound microscopic foundation yet. 

According to the spirit of the newly emerging economic discipline, a right 
direction to proceed is to reproduce observed data as they are. The analysis 
of financial data of listed firms in Japan shows that there are significant cor- 
relations among their capital, labor and value added; we regard those firms 
as constituting a statistical ensemble. Then we model such correlations in 
the production variables of firms using copulas (Sklar ( 1959[ )). Copulas are 
analytic functions created in statistics to describe genuine correlations among 
variables (Nelsen (2006), Joe (1997l). Mathematical properties of copulas 



have been studied extensively, although the concept of copulas has only a 
few decades of active history. Recently copulas have also attracted interest 
from business practitioners who have in their mind a number of possible ap- 
phcations to finance such as asset allocation, credit evaluation, default risk 



modeling, derivative pricing, and risk management (Cherubini et al. (2004|). 
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The copula model so obtained predicts value added yielded by a firm with 
given capital and labor in a probabilistic way. The new idea referred to as 
production copula is in sharp contrast to the production function which is 
completely deterministic. The present paper can thus be regarded as a gener- 
alization of the canonical economic concept at a microscopic level. And hence 
the production copula is expected to provide a tool to understand economic 
phenomena at both microscopic and macroscopic levels on an equal basis. 

In the following section we explain financial data on Japanese listed firms 
which are used in this paper. In Sec. Ill we fit the data to the production 
function of Cobb-Douglas form and demonstrate the real data cannot be ac- 
commodated in the framework of the production function. In Sec. IV we 
prepare the statistical modeling of production activities of firms in terms of 
copulas by briefly reviewing the concept of copulas. Section V is devoted 
to actual construction of the production copula. In Sec. VI we conclude this 
paper by pointing out possible applications of the production copula obtained. 



II. Data on Listed Firms in Japan 



The present analysis of real data is based on the NEEDS database (Nikkei 



Media Marketing, Inc. (2008)). It has accumulated financial statement data 
of listed firms in Japan over last 30 years; it is very exhaustive as regards 
financial information on the listed firms. We selected firms belonging to the 
manufacturing sector and compiled the fundamental quantities L, K, and 
Y for the production activities of those firms. To have a concrete idea, we 
actually substituted L, K, and Y respectively with labor cost, fixed asset, 
and value added; all of them are measured in units of million yen. Since 
labor cost L and fixed asset K are primary quantities, there is no ambiguity 
for them. However, there are basically two alternative ways (subtractive or 
additive method) to derive value added Y leading to different results. We 



refer the readers to Souma et al. ( 2009 ) for details about the calculation of Y 



We use the data set in 2006 throughout this paper; the total number of firms 
is TV = 13600 

In Fig.[T]we plot distributions for K , L, and Y of the manufacturers. These 
figures clarify the distribution of each quantity has a fat tail. We fitted these 
data adopting the following form, called the generalized beta function of the 
second kind (Kleiber and Kotz ( 2003| ), for the probability density function 
(PDF) p{x): 



p{x) 



B{tJ-/q,i^/q) X \xo^ 





1 + 




(xo) 











(2) 



where B(r, s) is the beta function. The corresponding cumulative distribution 
function (CDF) is expressed explicitly using the incomplete beta function 



^We have excluded one manufacturer with negative Y out of our database. 
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Figure 1: CDF's for the labor cost (a), the capital (b), and the 
value added (c) of the Japanese listed firms in the man- 
ufacturing sector, as depicted by dots. Results fitted to 
the real data in the form ^ are also shown by a solid 
line in each panel, along with the model parameters de- 
termined by the maximum likelihood method. 
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Figure 2: Correlations for each pair of the three financial quanti- 
ties, labor cost L, capital and value added Y , in the 
manufactures. 



B{z, r, s) as 



ith 



Py{x) = / dx'p{x') 



B{z,n/q,iy/q) 
B{fJ./q,i^/q) 



(3) 



(4) 



The parameters /i and v are power-law exponents in the large and small x 
limits, respectively: 





1 + 




z = 


(2^0) . 










for a; — > 00 , 



for X —f . 



(5) 



The parameter xg represents a characteristic scale and g is a crossover param- 
eter connecting the two limiting regions. Figure [1] confirms that Eq. ^ works 
well to reproduce the original data including the fat tail, where we determined 
the parameters using the maximum likelihood method. 

In Fig. [2] we show scatter plots for all pairs of the three financial quanti- 
ties in the manufactures. We see that those variables are mutually correlated, 
assuming that those firms constitute a statistical ensemble. Even L and K, 
which are regarded as independent variables in the production function, are 
correlated to a certain extent. A ridge theory of the production function elu- 
cidates that such dependence between L and K arises from profit maximiza- 



tion behaviors of firms (Aoyama et al. (2007), Souma (2007), Aoyama et al 



(2009)). Spearman's rank correlation coefficient pg is one of non-parametric 
measures of correlation between two variables. The values of ps for the three 
pairs, K-Y, L-Y, and L-K, are calculated as 0.86, 0.95, and 0.83, respectively. 
Thus the L-Y pair has stronger correlation than the remaining two pairs K-Y 
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and L-K have. And the remaining ones have correlations of similar strength. 
The detailed correlation structure for each pair is described in terms of copulas 
later, with an account of these results for the overall correlations. 



III. Production Function 

Before moving to the main theme in this paper, we spend short time here to 
discuss a production function appropriate to the data used here. One of the 
simplest functional forms for the production function ([T]) is given as 

F{L,K) = AK^'Lf^, (6) 

where A, a, and /3 are adjustable parameters to fit data. This form was intro- 



duced by Cobb and Douglasj ( 1928 1 a long time ago and has been extensively 
used because of its ease of use and its extreme flexibility. In addition we note 
that a generalized form, referred to as the CES production function, is also 



available (Arrow et al. (1961)). 

We have fltted the Cobb-Douglas (CD) form to the financial data of the 
listed firms using the least square method and have obtained the following 
values: 

A = 2.160, a = 0.183, /3 = 0.788 . (7) 

We note that the production function almost satisfies the homogeneity con- 
dition of degree one, namely, exhibits constant returns to scale: 

a + f3^ 0.971 . (8) 

The CD production function calculated with these parameters is compared 
with the original data in Fig. |3] Although one might be tempted to conclude 
that the CD production function describes the data well in viewing Fig. [3] it 
is far from the truth. 

Figure |4] shows difference between the actual value of Y and the value 
predicted by the best-fit CD production function F{L, K) with the parameters 
as given in Eq. ([8|. In the scatter plot (a) as well as in the histogram (b), 
we observe that the actual data is widely scattered around F{L,K). The 
standard deviation of this distribution is about 0.2, which means that the 
actual value of Y could easily be 10"-^ = 1.58, or 10^°'^ = 0.63 time the 
production function. The largest value of log{Y/ F{L, K)) is 1.78 {Y being 
60 times F{K,L)) and the smallest -1.09 (Y being 8% of F{L,K)). The 
existence of these extreme values implies the distribution of Y/F{K,L) has 
fat tails. 

Figure [5] shows this is in fact true; both for very small Y/F{L, K) and very 
large Y/ F{L, K). The cumulative distribution function (CDF) shows almost 
linear behaviors in log- log plot, which means that they obey power-laws. Be- 
cause of this existence of the fat tails on both sides, the production function 
often fails to describe the actual output Y by large ratio. While the produc- 
tion function describes only the central the central value of the distribution. 
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Figure 3: The CD production function (mesh plane) fitted to the 
real data (dots) in L-K-Y space viewed from two different 
angles. 



(a) (b) 




10^ K)-' 10* 10' 10« -1-0 -0.3 0.5 1.0 1 



F(L,K) \og(YIF(L,K)) 

Figure 4: Difference between the actual Y and the best-fit CD func- 
tion F{L,K); (a) scatter plot, where the dashed line 
shows the diagonal line, (b) histogram of number of firms 
in each bin of the logarithm of their ratio. 
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Figure 5: CDF of Y/F{L,K) in (a) [0, 1] region and (b) [l,oo] re- 
gion. 



it by no means can be used for any quantitative theory involving Y . Instead 
we need a theoretical device that allows us to describe the distribution itself, 
which we will consider in the following sections. 

IV. Copula Theory 

A natural way to incorporate the stochastic nature of Y even at given K and 
L is to construct a probability density function (PDF) model for the three 
financial quantities. The PDF is defined as 

1 ^ 

p(L,i^,r) = -^5(L-L,),5(i^-i^,),5(r-r,), (9) 

1=1 

where the quantities with a subscript i are those of a firm i and 6{x) denotes 
Dirac's delta function. The CDF corresponding to p {L, K,Y) is then given 

by 

pL nK nY 

P<{L,K,Y)= / / / dL'dK'dY' p{L' ,K' ,Y'). (10) 
Jo Jo Jq 

The marginal PDF's with two variables are deduced from p{L, K, Y) by inte- 
grating over one of the variables, i.e., 

p{L,Y)= f dKp{L,K,Y) . (11) 

This relation is cast onto the CDF's as 

P<{L,Y)= f f dL'dY'p{L',Y') = P<{L,K^oo,Y) . (12) 
Jo Jo 

Then the marginal PDF's for each variable are likewise deduced from the 
binary PDF's such as 



p{Y) 



[ dLp{L,Y), (13) 
Jo 
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for instance. The marginal CDF is thus related to the binary CDF through 



P<(y) = / dY'p{Y') = P<{L^ oo, Y) . (14) 
Jo 

The PDF of Y at given L and K is given as a conditional probability: 

piY\L,K)^P-^^. (15) 

This is an exact description of the data: unlike the production function ([T]), 
Y is not determined uniquely as a function of K and L, but obeys a PDF 
whose parameters are determined by the given K and L. Our task is thereby 
to adopt a mathematical measure that is suitable for describing this PDF, 
which we will do next. 



A. Definition of copula 



To make an analytic model for the PDF, p{L, K,Y), we take advantage of 
a copula method. Copula is statistical means to measure dependence among 
stochastic variables. In other words, copulas extract correlations inherent 
among stochastic variables free from the marginal CDF's of variables them- 
selves. A number of workable forms for copulas have been proposed by statis- 
ticians along with elucidation of their mathematical properties. 



Sklar's theorem ( |Sklar| ( |1959| ); |Nelsen| ( |2006| )) guarantees that the CDF 
P^{L, K,Y) is a unique function of the marginal CDF's associated with L, 
K, and Y: 

P<mK,Y)^C{uL,UK,UY), (16) 

where 

Us -.^ P<{s) = 1 - Py{s) {s = L,K,Y). (17) 

The function C{ul,uk ,uy) is called copula. The PDF, p{L,K,Y) is then 
derived from P^(L, K,Y) by carrying out partial differentiation with respect 
to each of the variables: 



p{L,K,Y) 



d^P<{L,K, Y) 



p{L)p{K)p{Y)c {ul.uk, uy) 



dLdKdY 

where we introduced the copula density c(u/,, u^, uy) defined by 



c{ul,uk,uy) 



d^C {ul,uk,uy) 
duLduxduY 



(18) 



(19) 



If the variables are statistically independent of each other, then the copula 
density reduces to 1. Thus the copula density, referred to as the correlation 



function in many-body physics ( |Hansen and McDonald (2006 1), and hence 
the copula describe genuine correlations among the variables. 
The copulas have the boundary conditions exemplified as 



C {ul,uk = 0,-uy) = 



(20) 
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C{ul,uy) = C {ul,uk = l,wy) 



(21) 



The condition ( 20 1 is transparent from the definition ( 10 1 for the CDF; Eq. ( 21 1 
just rephrases the relation (12) (Nelsen ( 2006[ )). Also we note that if the fi- 
nancial quantities are independent of each other, then the copulas read, 



C{Us,Us') = UsUs', 

C{ul,uk,uy) = ulUkUy- 
B. Archimedean copulas 



(22) 
(23) 



There is a well-known family of copulas called Archimedean copulas, which 
have been widely used because of their ease of mathematical handling together 
with diversity of their correlation properties. We first discuss the Archimedean 
copulas with two variables and then proceed to those with many variables. 

A bivariate Archimedean copula can be readily constructed with a gener- 
ator function rj (u) as 



Ca{ui,U2) = T] ^ (77(1*1) + 77(^2)) . 



(24) 



Here we assume that 'ri{u) is a continuous and strictly decreasing function 
mapping [0,1] to [0,oo] such that 77(1) = 0. We also have to impose convexity 



on ri(u) so that the constructed function (24) is a valid copula 



Here we test three typical Archimedean copulas to model the real data. 
One of them is the Frank copula given by 



Vf{u;9) 
Cf {ui,U2;6) 



-Bu 



log 



4 log 



1 



1 



- 1 

(e-8"i 



with 
-l)(e 



' CO < < 00, 
-eu2 _ i\ 



(e-«-l) 

The second one is referred to as the Gumbel copula: 

■q q{u; 9) — {—log uf with 1 < 6* < 00, 



Cg {ui,U2\ 



exp 



|(-log7il)'' + (-logTiz)*'! 



1/6' 



The last one is the Clayton copula given by 
77c(u;^)= ^(t.-'^-I) 
Cc{ui,U2\0) — [max (u'^^ + -u 



with - 1 < < cx), 
-i/e 



i>o)]- 



(25) 
(26) 

(27) 

(28) 

(29) 
(30) 



These three copulas possess correlation properties of different character- 
istics as displayed in Fig. |6] The correlation strength of the Frank copula 
is almost flat spanning from the lower tail (ui ,U2 — 0) to the upper tail 
(ui,M2 ~ 1). The Gumbel copula has stronger correlation in the upper tail 
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Figure 6: Contour plots of the binary copula densities. The panel 
(a) is for Frank with 9 = 14.14; (b), for Gumbel with 
= 4; (c), for Clayton with 6 = 6. The parameter 6 is 
chosen to yield the same value of the Kendall t (= 0.75) 
for each copula. 



than in the lower tail and the Clayton copula shows reversed correlation struc- 
ture. This difference in asymmetry in correlation structure is, however, su- 
perficial. In fact, by changing the variables, Ui — > 1 — u^, we derive a different 
copula model from a given copula, called survival copula. The survival C lav- 
ton copula (s-Clayton, for short), Cq {ui,U2]9), is related (Nelsen (2006 1) to 
the original Clayton copula through 



Cc (ui, U2; 6*) = ui + U2 - 1 + Cc (1 - ui, 1 - U2; 



(31) 



As is easily appreciated, in turn, Cc(ui,U2) has stronger correlation in the 
upper tail. 

The bivariate Archimedean copulas can be generalized to multivariate cop- 
ulas with n variables in such an iterative way as 



Ca (ui, • • • ,u„) = 77 ^(77(mi)H h77(u„)) 



(32) 



The generator function however is required to satisfy a more stringent math- 
ematical constraint that its inverse ri~^{x) is completely monotonic, i.e.. 



(-1)' 



>0 (m = 0, 1,2,---). 



(33) 



We refer the readers to the textbook by NelsenJ ( 2006 ) for the details. The 



generators of Frank, Gumbel, and Clayton are completely monotonic with 
the following conditions for 9: 9 > (Frank), 9 > 1 (Gumbel), and 9 > 
(Clayton). We will call a multivariate Archimedean copula after its generator; 



namely, the copula ( 32 1 generated by Eq. ( |25[ ) is referred to as n-variate Frank 
copula. 



C. Generalization of Archimedean copulas 

The Archimedean copulas, widely used for applications, may be sometimes 
too restrictive to accommodate real data. 
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They are totally symmetric with respect to exchange of independent vari- 
ables. A way to relax this constraint for the bivariate copulas is to set 



'•'2 ; ' 



(34) 



where two additional parameters with < a, /? < 1 are introduced. In the 
case oi a = f3 the generalized copula is exchangeable. Difference between 
a and f3 thus measures the degree of asymmetry in correlations. Taking a 
special limit of a ~ /3 — 1 recovers the original Archimedean copulas. 



Also we recall that the trivariate Archimedean copula, Eq. ( 32 1 with n — 3, 
is characterized by a single parameter 6. It means that all of the marginal 
CDF's have the same correlation structure. This is not true for the real data 
under study as has been already indicated, so that we generalize the original 
form as 



CA-ncK{ui,U2,U3;9i,62) = Ca (Ca (ui , U2 ; ^2 ) , W3 ; ^'l ) 



= Vi ^ {vi {v2 ^ iv2 (ui) + m (u2))) + m (us)) 



(35) 



with different generator functions rji^u) = r]{u;9i) and r]2{u) = rj{u;92)- 



The copula ( 35 1 now has two characteristic parameters 9i and 92 leading 



to the margin Ca(ui, 1*2; ^i) generated by 77i(m) and the remaining margins, 
Ca(ui, U3; 6*2) and Ca(u2, M3; 6*2) generated by ri2{u). If one set 9i = 02 in Eq. 
( 35 1 , the exchangeable Archimedean copula ([32| with n = 3 is recovered. We 



note that the condition 9i > 92 should be satisfied for Eq. ( 35 1 to be a cop- 



ula. We refer to the generalized form (351 as non-exchangeable Archimedean 



copula following McNeil et al. (20051 



V. Construction of Production Copula 

In this section the distributed real data are modeled using copulas as they 
are. We refer to the resulted copula as "production copula" . 



A. Bivariate 

First we pay our attention to binary correlations for all pairs out of L, K, 
and Y as are depicted in Fig. |2] Three Archimedean copulas due to Frank, 
Gumbel and Clayton are adopted. 

We have carried out the maximum likelihood estimate to fit those copulas 
to the pair correlations as given in Fig. [2] The fitting results are summarized 
in Table [l] where the optimized values of the correlation parameter and 
the maximized log-likelihood I are listed. Also we calculated the Kendall 
rank correlation coefficient r from the 6''s for each copula. We observe that 
the Gumbel copula gives the best fitting among the three copulas tested as 
manifested by the largest value for The optimized results with the Gumbel 
copula are thus compared with the corresponding empirical copulas derived 
from the real data in Fig. |7] The two results are in good agreement with each 
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Table 1: Maximum likelihood estimate of bivariate Archimedean 
copulas fitted to the pair correlations in Fig. |2] Also the 
Kendall rank correlation coefficient r is listed correspond- 
ing to the value of 9 for each copula. 



Copula 


K-Y 


L-Y 


L-K 


Frank 

i 

T 


11.0 
940.1 
0.691 


21.2 
1612.4 
0.826 


9.23 
790.2 
0.644 


Gumbel d 

i 

T 


3.21 
1081.8 
0.688 


5.30 
1694.0 
0.811 


2.73 
892.1 
0.634 


s-Clayton d 

I 

T 


3.43 
992.1 
0.632 


6.11 
1483.7 
0.753 


2.59 
787.3 
0.564 




Figure 7: Modeling of the pair correlations in Fig.[2|in terms of the 
Gumbel copula. The fitted results (smooth solid curves) 
are compared with the corresponding real data (jaggy 
solid curves) ; difference between the two curves is almost 
invisible. 



other. The same comparison but for the copula densities are also made in 
Fig.|8] 

The overall strength of the pair correlations is represented by magnitude of 
9 and r. Relative comparison of 9's and t's in Table [T] confirms the conclusion 
drawn in the previous section. The L-Y correlation is significantly stronger 
than the K-Y and L-K correlations. In contrast the last two correlations 
resemble each other. 

In passing we address to what extent the binary correlations among L, K, 
and Y are asymmetric. We repeated the maximum likelihood estimate using 
the asymmetric Gumbel copula given by Eq. (|34| with Eq. ([27]). Tabled 
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Figure 8: Contour plots of the copula densities corresponding to 
the fitted copulas in Fig. |7] The dots refer to the real 
data. 



Table 2: Degree of asymmetry in the pair correlations in Fig. [2| 
estimated using the asymmetric Gumbel copula. 





K-Y 


L-Y 


L-K 


6 


3.34 


6.59 


2.86 


a 


0.994 


0.999 


1.000 




0.984 


0.934 


0.957 


i 


1090.1 


1779.5 


896.0 



demonstrates that none of the pairs has notable asymmetry in its correlation 
structure. 

B. Trivariate 

At last we are ready to discuss the principal objective of the present paper, 
namely, construction of the production copula. The models we consider here 
are enumerated as 

Model (I) Frank copula; 

C'^^\ul,uk,uy) = Cf{ul,uk,uy;0). (36) 
Model (II) Gumbel copula; 

C'^''Hul,uk,uy) = CGiuL,UK,UY;e). (37) 

Model (III) s-Clayton copula; 

C'^'''\uL.UK,UY) = CciuL ,UK,UY\a) 

= ul + uk + uy ~ I + Cc(l - ul, 1 — uk^ 1 — uy;0). 

(38) 
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Model (IV) Non-exchangeable Gumbel copula; 

{ul,uy;02),uk;0i). (39) 



In the above, trivariate Archimedean copulas Cp^Q^c{uL,uji,UY',0) are con- 
structed from their 77-functions ( 25 1 , ( 27 ) , ( 29 ) respectively by the prescription 
( 32 1 . The first three models specified by a single parameter result in all of 
marginal CDF's and PDF's with the same correlation structure. The last 
model is based on Eq. (35 1 having Eq. (28 1 plugged into Ca (ui , ; 6*1 ) and 
CA(ui, Uj-; 6*2). The variables ui, M2, and 1*3 in Eq. (351 should read u^, uy, 
and uk, respectively. The last model can account for the difference in cor- 
relations among the marginals with different values for 9i and 62, which is a 
desirable property as we will see. 



Table 3: Maximum likelihood optimization in the models for the 
production copula. 





Model (I) 


Model (II) 


Model (III) 


Model (IV) 


e 


11.5 


3.27 


3.34 


2.89 (0i), 5.26 {02) 


£ 


2207.5 


2428.2 


2164.3 


2701.1 


AIC 


-4413.0 


-4854.4 


-4326.6 


-5398.2 



We adapted the four models to the real data using the maximum likeli- 
hood method again. Results of the optimization are listed in Table [3] As has 
been easily expected from the experience in modeling the binary correlations, 
the Gumbel copula used in Model (II) certainly outperforms the fitting as 
compared with the other Archimedean copulas in Models (I) and (III). This 
is actually the motivation behind Model (IV) in which the Gumbel copula 
is specially selected. The result of the optimization based on the general- 
ized model is also included in Table [3] Since the numbers of parameters are 
different between Models (II) and (IV) , we have to replace the maximum log- 



likelihood £ by Akaike's information criterion AIC for model selection ( Akaike 



(1974)); a model with the smallest AIC should be adopted. The criterion is 
given by 

AIC = -2£ + 2k, (40) 

where k is the number of parameters in a statistical model. Comparison of 
the AIC values proved that Model (IV) significantly improves the fitting over 
even Model (II). 

Further we delve into the performance of Model (IV). To make detailed 
comparison with the real data, we employ a trivariate copula cumulani|^ de- 

^Cumulants are usually defined in terms of the PDF's in place of the CDF's ( [Hansen] 
[and McDonald| l |2006| l). 
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fined by 



n{uL,UK, Uy) = C{ul,Uk, Uy) - UkC{ul,Uy) 

— uyC{ul,uk) + 2ulukUy- 



UlC{uk,Uy) 



(41) 



Here contributions essentially due to binary correlations are subtracted from 
the trivariate copula. To appreciate this fact, we will take two special cases. 
First suppose that all of the variables are independent of each other. Replace- 



ment of the copulas on the right-hand side of Eq. ( 41 1 by the corresponding 
independent copulas (22 1 and (23 1 leads to n{uL,UK,UY) = 0. Next suppose 
that only a pair of L and K are intercorrelated among the three variables, 
then the copula C{ul, uk, uy) is decomposed as 



C{ul,uk,uy) = uyC{ul,uk)- 



(42) 



And the bivariate copulas involving uy as a variable can be replaced by the 
corresponding independent copulas. Again we see that ^{uj^, uk, uy) vanishes 
in this case. Such subtraction is also manifested by the boundary condition 
that the copula cumulant vanishes on the marginal boundaries, namely. 



n{l,UK, Uy) = fi(uL, 1, Uy) = ^{ul,uk, 1) = 0. 



(43) 



Figure |9] shows the fitted results for ^{ul,uk,uy) in Model (IV) together 
with the corresponding empirical data on the typical cross sections A-D de- 
picted in Fig. [TOj We observe the copula model reproduces the empirical 
results in an almost indistinguishable manner. Figure [TT] gives more in-depth 
comparison of the results between the copula model and the real data for 
U,{x,x,x) along the diagonal direction specified hy x ~ u^ = uk = uy. 
Again switching from Model (II) to (IV) leads to significant improvement in 
reproducing the empirical copula cumulant. 

Now that the copula Model (IV) has been established, it may be used for 
various studies of firms' activities. Here we demonstrate one such an example. 
We carried out a simulation of production activities of the Japanese listed 



firms in 2006 and exhibit the result in Figure 12 The simulated points {N 



1360) were generated using the rejection method (Press et al. (2007 1 ) with the 
trivariate PDF obtained by combining the production copula and the marginal 
PDF's according to Eq. (18 1. At the same time we did a simulation based on 



the random model where no correlations among the financial quantities were 



taken into account; the result is added in Fig. 12 Those results are compared 
with the real data depicted in Fig. |3] 

Also we calculated the deviation from the CD production function for the 
copula and random models, as has been already done for the real data in 
Figs. |4]and[5] One can express the two CDF's of the ratio ^ = Y / F{L, K) 
in terms of the trivariate copula. For instance, the CDF in the upper side is 
given as 

P>{0^f{C}/f{l) (l<^<oo), (44) 
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Real data 



Copula model 




UK 



UK 




D 



1 

Uy 






ul 1 




■ 







Ul 



1 




Figure 9: Contour plots of the copula cumulant, Eq. (41 1 on the 
cross sections in u^-uk-uy space such as depicted in 
Fig. [lOj The results obtained from the real data (left- 
hand side) are compared with those derived from the 
copula model (right-hand side); the contours are drawn 
at the same levels on both sides. 
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1 




Figure 10: Typical cross sections in ui^-uk-uy space. 
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Figure 11: Comparison of the results for Q.{x, x, x) between the 
copula model (solid curve) and the real data (dots). 
Model (II) is adopted in the panel (a); Model (IV), in 
the panel (b). 
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(a) copula model 




(b) random model 




Figure 12: A simulated result for production activities of the 
Japanese listed firms obtained with the production cop- 
ula in Model (IV), accompanied by that based on the 
random model without any correlations among the fi- 
nancial quantities. The 3D graphs are drawn with the 
same view angles as in Fig.lSlfor each model. 
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Figure 13: CDF of Y/F{L, K) calculated with the copula Model 
(IV) compared to the results based on the real data 
and the random model, corresponding to Fig. Isj 



where 

d'^C{uL,UK,UY,) 



fiO = 1 - / dUL [ du 

Jo Jo 



^ du^duK 



(45) 



with 

= iF{L,K) = CF{P-\ul),P^\uk)) • (46) 



For the random model, Eq. ( 45 ) is replaced by 



/(O = 1 - / duL / duKUY^ ■ 
Ja Jo 



(47) 



Detailed comparison of those results with the corresponding real data is made 
in Fig. 13 The analytic model ([3| for the marginal CDF's was highly useful 



in executing the integration in Eq. ( |45| . We see that correlations involved 
in the real data are taken into account accurately by the copula model, ex- 
cept for the outliers in the upper side occupying about 5% of the data; the 
functional behavior of the outliers is rather close to that of the results in the 
random model. We thus claim that the production copula is so successful in 
reproducing the real data and so workable considering that it has only two 
adjustable parameters. 



VI. Conclusion 



In this paper we have proposed the use of the production copula to take full 
account of a wide variety of production activities of firms. 

We showed that capital, labor and value added of firms are closely related 
in a probabilistic way to each other by analyzing financial data of listed firms 
in Japan. At the same time we confirmed that the productive heterogeneity 
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of firms was far beyond the scope of a production function. These empirical 
facts authorizes our endeavor to model the production activities of firms in 
terms of copulas. Four copula models were employed and their accuracy was 
examined through fitting to the real data. The production copula so obtained 
predicts the value added yielded by a firm with given capital and labor in a 
probabilistic way. This gives rise to paradigm shift in the economic theory 
for studying production activities of firms; generalization of such a canonical 
economic concept as the production function at a microscopic level. 

To demonstrate possible applications of the production copula, we carried 
out a simulation for production activities of the Japanese listed firms. Also we 
noted that the production copula would enable us to predict the GDP with 
statistical uncertainty taking into account the diversity in firms' productivity. 

We thus believe that the production copula will play a role of Jacob's 
Ladder to go back and forth between micro and macro economics. With 
the production copula, for instance, it is not necessary to solve such a long- 
standing aggregation problem as microscopic foundation of the macroscopic 
production function. 
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